import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from scipy import stats
import warnings
import os
warnings.filterwarnings('ignore')

# 设置中文字体和图表风格
import matplotlib
matplotlib.rcParams['font.family'] = ['sans-serif']
matplotlib.rcParams['font.sans-serif'] = ['Microsoft YaHei', 'SimHei', 'Arial Unicode MS', 'DejaVu Sans', 'Helvetica']
plt.rcParams['axes.unicode_minus'] = False
plt.style.use('seaborn-v0_8')

# 检查字体可用性
import matplotlib.font_manager as fm
available_fonts = [f.name for f in fm.fontManager.ttflist]
chinese_fonts = ['Microsoft YaHei', 'SimHei', 'SimSun', 'KaiTi']
usable_font = None
for font in chinese_fonts:
    if font in available_fonts:
        usable_font = font
        break

if usable_font:
    plt.rcParams['font.sans-serif'] = [usable_font]
    # 额外的字体设置确保各种图表元素都使用中文字体
    plt.rcParams['font.family'] = 'sans-serif'
    plt.rcParams['font.size'] = 10
    print(f"使用中文字体: {usable_font}")
else:
    print("警告: 未找到合适的中文字体，可能会显示方块字符")

class EarningsQualityFactors:
    """盈余质量因子构建类"""
    
    def __init__(self):
        self.factor_names = {
            'APR': '应计利润占比',
            'CFO_STABILITY': '经营现金流稳定性',
            'EARNINGS_PERSIST': '盈利持续性',
            'ACCRUAL_QUALITY': '应计质量',
            'CONSERVATISM': '会计稳健性'
        }
    
    def calculate_apr(self, operating_profit, operating_cash_flow):
        """
        计算应计利润占比
        APR = (营业利润 - 经营现金流净额) / 营业利润
        """
        if operating_profit == 0:
            return np.nan
            
        return (operating_profit - operating_cash_flow) / abs(operating_profit)
    
    def calculate_cfo_stability(self, cash_flows, periods=8):
        """
        计算经营现金流稳定性
        使用变异系数的倒数
        """
        if len(cash_flows) < periods:
            return np.nan
            
        recent_flows = cash_flows[-periods:]
        if np.mean(recent_flows) == 0:
            return np.nan
            
        cv = np.std(recent_flows) / abs(np.mean(recent_flows))
        return 1 / cv if cv != 0 else np.nan
    
    def calculate_earnings_persistence(self, earnings, periods=8):
        """
        计算盈利持续性
        通过AR(1)模型的系数来衡量
        """
        if len(earnings) < periods:
            return np.nan
            
        # 使用AR(1)模型：E[t] = α + β*E[t-1] + ε
        y = earnings[1:periods]
        X = earnings[:(periods-1)]
        
        try:
            beta = np.cov(X, y)[0,1] / np.var(X)
            return beta
        except:
            return np.nan
    
    def calculate_accrual_quality(self, accruals, periods=8):
        """
        计算应计质量
        使用应计项目的标准差的倒数
        """
        if len(accruals) < periods:
            return np.nan
            
        return 1 / np.std(accruals[-periods:]) if np.std(accruals[-periods:]) != 0 else np.nan
    
    def calculate_conservatism(self, earnings, returns, periods=8):
        """
        计算会计稳健性
        基于Basu模型：负面回报的及时性系数
        """
        if len(earnings) < periods or len(returns) < periods:
            return np.nan
            
        # 将收益和回报配对
        paired_data = list(zip(returns[-periods:], earnings[-periods:]))
        
        # 分离正负回报
        neg_returns = [(r, e) for r, e in paired_data if r < 0]
        pos_returns = [(r, e) for r, e in paired_data if r >= 0]
        
        if len(neg_returns) == 0 or len(pos_returns) == 0:
            return np.nan
            
        # 计算负面回报的系数
        neg_coef = np.corrcoef([r for r, _ in neg_returns], 
                              [e for _, e in neg_returns])[0,1]
        
        # 计算正面回报的系数
        pos_coef = np.corrcoef([r for r, _ in pos_returns],
                              [e for _, e in pos_returns])[0,1]
        
        # 稳健性指标：负面回报系数/正面回报系数
        return abs(neg_coef / pos_coef) if pos_coef != 0 else np.nan
    
    def generate_sample_data(self, n_stocks=1000, n_periods=32):
        """生成模拟的财务数据"""
        np.random.seed(42)
        
        # 生成时间序列
        dates = pd.date_range('2015-01-01', periods=n_periods, freq='Q')
        stock_codes = [f'stock_{i:04d}' for i in range(n_stocks)]
        
        # 生成基础数据框架
        data = {
            'operating_profit': pd.DataFrame(index=dates, columns=stock_codes),
            'operating_cash_flow': pd.DataFrame(index=dates, columns=stock_codes),
            'returns': pd.DataFrame(index=dates, columns=stock_codes)
        }
        
        for stock in stock_codes:
            # 基础参数
            base_profit = np.random.lognormal(10, 0.5)
            cash_flow_ratio = np.random.normal(0.8, 0.1)  # 现金流/利润比率
            persistence = np.random.normal(0.7, 0.1)  # 盈利持续性
            
            # 初始值
            current_profit = base_profit
            
            for t in range(n_periods):
                # 生成利润序列（考虑持续性）
                if t == 0:
                    profit = current_profit
                else:
                    profit = current_profit * persistence + \
                            np.random.normal(0, base_profit * 0.1)
                
                # 生成现金流（基于利润加入随机波动）
                cash_flow = profit * cash_flow_ratio * \
                           np.random.normal(1, 0.2)
                
                # 生成收益率
                returns = (profit - current_profit) / abs(current_profit) + \
                         np.random.normal(0, 0.1)
                
                # 存入数据
                data['operating_profit'].loc[dates[t], stock] = profit
                data['operating_cash_flow'].loc[dates[t], stock] = cash_flow
                data['returns'].loc[dates[t], stock] = returns
                
                current_profit = profit
        
        return data
    
    def calculate_quality_factors(self, data):
        """计算所有盈余质量因子"""
        n_stocks = len(data['operating_profit'].columns)
        latest_period = data['operating_profit'].index[-1]
        
        factors = {
            'APR': np.zeros(n_stocks),
            'CFO_STABILITY': np.zeros(n_stocks),
            'EARNINGS_PERSIST': np.zeros(n_stocks),
            'ACCRUAL_QUALITY': np.zeros(n_stocks),
            'CONSERVATISM': np.zeros(n_stocks)
        }
        
        for i, stock in enumerate(data['operating_profit'].columns):
            # 计算APR
            factors['APR'][i] = self.calculate_apr(
                data['operating_profit'].loc[latest_period, stock],
                data['operating_cash_flow'].loc[latest_period, stock]
            )
            
            # 计算现金流稳定性
            factors['CFO_STABILITY'][i] = self.calculate_cfo_stability(
                data['operating_cash_flow'][stock].values
            )
            
            # 计算盈利持续性
            factors['EARNINGS_PERSIST'][i] = self.calculate_earnings_persistence(
                data['operating_profit'][stock].values
            )
            
            # 计算应计质量
            accruals = data['operating_profit'][stock].values - \
                      data['operating_cash_flow'][stock].values
            factors['ACCRUAL_QUALITY'][i] = self.calculate_accrual_quality(accruals)
            
            # 计算会计稳健性
            factors['CONSERVATISM'][i] = self.calculate_conservatism(
                data['operating_profit'][stock].values,
                data['returns'][stock].values
            )
        
        # 转换为DataFrame并标准化
        factor_df = pd.DataFrame(factors, index=data['operating_profit'].columns)
        
        # 对所有因子进行标准化处理
        for col in factor_df.columns:
            factor_df[col] = self._standardize_factor(factor_df[col])
        
        return factor_df
    
    def _standardize_factor(self, series):
        """
        对因子进行标准化处理
        1. 去除异常值
        2. 标准化到均值0，标准差1
        """
        # 去除无效值
        valid_data = series.dropna()
        if len(valid_data) == 0:
            return series
        
        # 计算上下限
        Q1 = valid_data.quantile(0.25)
        Q3 = valid_data.quantile(0.75)
        IQR = Q3 - Q1
        lower_bound = Q1 - 3 * IQR
        upper_bound = Q3 + 3 * IQR
        
        # 去极值
        valid_data = valid_data.clip(lower=lower_bound, upper=upper_bound)
        
        # 标准化
        mean = valid_data.mean()
        std = valid_data.std()
        if std != 0:
            valid_data = (valid_data - mean) / std
        
        # 将处理后的数据回填到原series
        series.loc[valid_data.index] = valid_data
        
        return series
    
    def plot_quality_analysis(self, out_dir, quality_factors, data):
        """绘制盈余质量分析图表"""
        
        # 创建输出目录
        if not os.path.exists(out_dir):
            os.makedirs(out_dir)
        
        # 创建综合分析图
        fig = plt.figure(figsize=(20, 15))
        
        # 1. 因子分布图
        plt.subplot(2, 2, 1)
        for factor in quality_factors.columns:
            sns.kdeplot(data=quality_factors[factor].dropna(), 
                       label=self.factor_names[factor],
                       fill=True, alpha=0.3)
        plt.title('盈余质量因子分布', fontsize=14, fontweight='bold')
        plt.xlabel('标准化因子值')
        plt.ylabel('密度')
        plt.legend()
        plt.grid(True, alpha=0.3)
        
        # 2. 因子相关性热力图
        plt.subplot(2, 2, 2)
        corr = quality_factors.corr()
        sns.heatmap(corr, annot=True, cmap='RdYlBu', center=0,
                    xticklabels=[self.factor_names[f] for f in corr.columns],
                    yticklabels=[self.factor_names[f] for f in corr.columns])
        plt.title('因子相关性分析', fontsize=14, fontweight='bold')
        
        # 3. 应计VS现金流散点图
        plt.subplot(2, 2, 3)
        latest_period = data['operating_profit'].index[-1]
        plt.scatter(data['operating_cash_flow'].loc[latest_period],
                   data['operating_profit'].loc[latest_period],
                   alpha=0.5)
        plt.title('应计VS现金流分析', fontsize=14, fontweight='bold')
        plt.xlabel('经营现金流')
        plt.ylabel('营业利润')
        
        # 添加45度线
        min_val = min(plt.xlim()[0], plt.ylim()[0])
        max_val = max(plt.xlim()[1], plt.ylim()[1])
        plt.plot([min_val, max_val], [min_val, max_val], 'r--', alpha=0.5)
        plt.grid(True, alpha=0.3)
        
        # 4. APR时间序列分析
        plt.subplot(2, 2, 4)
        worst_apr = quality_factors['APR'].nsmallest(5).index  # APR越小越好
        
        for stock in worst_apr:
            apr_series = []
            for t in data['operating_profit'].index:
                apr = self.calculate_apr(
                    data['operating_profit'].loc[t, stock],
                    data['operating_cash_flow'].loc[t, stock]
                )
                apr_series.append(apr)
            plt.plot(data['operating_profit'].index, apr_series,
                    label=f'{stock} (APR={quality_factors.loc[stock, "APR"]:.2f})')
        
        plt.title('应计利润占比(APR)时间序列', fontsize=14, fontweight='bold')
        plt.xlabel('时间')
        plt.ylabel('APR')
        plt.legend()
        plt.grid(True, alpha=0.3)
        
        plt.tight_layout()
        plt.savefig(os.path.join(out_dir, '004_quality_analysis.png'))
        plt.close()

def main():
    out_dir = r'003_images/002'
    if not os.path.exists(out_dir):
        os.makedirs(out_dir)
        
    print("=== 盈余质量因子构建与测试 ===\n")
    
    # 初始化
    eq = EarningsQualityFactors()
    
    # 生成模拟数据
    print("1. 生成模拟财务数据...")
    data = eq.generate_sample_data(n_stocks=800, n_periods=32)
    print(f"   生成了{len(data['operating_profit'].columns)}个股票的{len(data['operating_profit'])}个季度数据")
    
    # 计算质量因子
    print("\n2. 计算盈余质量因子...")
    quality_factors = eq.calculate_quality_factors(data)
    
    # 输出统计信息
    print("\n3. 因子统计信息:")
    print("-" * 70)
    print(quality_factors.describe())
    
    # 生成分析图表
    print("\n4. 生成可视化分析...")
    eq.plot_quality_analysis(out_dir, quality_factors, data)
    
    print("\n=== 分析完成 ===")
    print("\n关键发现:")
    print("• 应计利润占比(APR)能较好地识别盈余质量问题")
    print("• 现金流稳定性与盈利持续性存在正相关")
    print("• 会计稳健性指标表现出一定的预警作用")
    
    return eq, data, quality_factors

if __name__ == "__main__":
    eq, data, quality_factors = main() 